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Abstract. The effects of elasticity on the break-up of liquid threads in microfluidic cross-junctions is investigated using 
numerical simulations based on the "lattice Boltzmann models" (LBM). Working at small Capillary numbers, we investigate 
the effects of non-Newtonian phases in the transition from droplet formation at the cross-junction (DCJ) and droplet formation 
downstream of the cross-junction (DC) (Liu & Zhang, Phys. Fluids. 23, 082101 (2011)). Viscoelasticity is found to influence 
the break-up point of the threads, which moves closer to the cross-junction and stabilizes. This is attributed to an increase of 
the polymer feedback stress forming in the corner flows, where the side channels of the device meet the main channel. 
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INTRODUCTION 

The progressive breakup of a fluid thread into a number of small drops is a rich physical phenomenon [1]. It is also 
an important dynamical process that impacts many applications. In particular, droplet-based microfluidic devices have 
gained a considerable deal of attention due to their importance in studies that require high throughput control over 
droplet size [2, 3, 4, 5]. Common droplet generator designs used in these devices are T-shaped geometries [6] and 
flow-focusing devices [7, 8, 9]. In T-shaped geometries, a dispersed (d) phase is injected perpendicularly into the main 
channel containing a continuous (c) phase. Forces are created by the cross-flowing continuous phase which periodically 
break off droplets; the flow-focusing geometry, instead, creates drops by focusing a fluid stream into another co¬ 
flowing immiscible fluid [10, 11, 12, 7, 8]. The operational regime of these devices is primarily characterized by the 
Capillary number Ca, which quantifles the importance of the viscous forces with respect to the surface tension forces 
at the non-ideal interface and the drop size and its distribution are dictated by the relative flow rates of the two fluids. 
With few exceptions [13, 14, 15], previous research has been mainly restricted to Newtonian fluids. The processing 
of biological fluids with high molecular weight macromolecules, inevitably results in considering a non-Newtonian 
(NN) behaviour. Consistently, the use of viscoelastic liquids in flow-focusing devices [13, 16, 14] or T-junction 
geometries [15] has recently gained attention. The formation and the pinch-off mechanism of viscoelastic droplets 
in Newtonian continuous phases was investigated in various geometries of the flow-focusing device by Steinhaus et 
al [14]. The effect of the molecular weight of the polymer in the droplet phase on filament thinning was studied by 
Arratia et al. [13, 16]. Typically, it is observed that elasticity of the droplet liquid prolongs the processes of thinning 
of the filament and significantly increases the interval required for breakup to complete. In a recent paper Garstecki 
et al [9] presented an experimental comparative study of the effect of elasticity of the continuous liquid in generation 
of droplets in microfluidic flow-focusing devices. The authors And that the elasticity of the focusing liquid stabilizes 
the jets facilitating formation of smaller droplets, and leads to transitions between various regimes at lower ratios of 
flow and at lower values of the Capillary numbers in comparison to the Newtonian focusing liquids. Complementing 
these results with systematic investigations by varying deformation rates and fluid constitutive parameters would be of 
extreme interest. This is witnessed by the various papers in the literature [6, 7, 8] addressing these kind of problems 
with the help of numerical simulations. In a recent paper, Liu & Zhang [7,8] performed lattice Boltzmann simulations 
of three-dimensional microfluidic cross-junctions at low Capillary numbers. A regime map was created to describe 
the transition from droplets formation at the cross-junction (DCJ), downstream of the cross-junction (DC), to stable 
parallel flows (PF). The influence of flow rate ratio. Capillary number, and channel geometry was then systematically 
studied in the squeezing-pressure-dominated DCJ regime. In the present paper we study the impact of NN continuous 
phases in the transition from DCJ to DC regime. 


THEORETICAL MODEL 


Our numerical approach is based on a hybrid combination of LBM and finite difference schemes, the former used to 
model two immiscible fiuids with a variable viscous ratio, and the latter used to model viscoelasticity using the FENE-P 
constitutive equations [26, 27]. LBM have already been used to model droplet deformation problems [17, 18, 19, 20], 
droplets in confined T-junctions [20, 21] and also viscoelastic fiows [22, 23]. This approach has been studied and 
validated in other dedicated works [24, 25], where we have provided evidence that the model is able to capture 
quantitatively rheological properties of dilute suspensions as well as deformation and orientation of single droplets 
in confined shear fiows. We just recall here the main reference equations which are integrated in both the continuous 
and dispersed phases. 

In the continuous phase we integrate both the NS (Navier-Stokes) equations for the velocity Uc and the EENE-P 
equations 


pc[dtUc + {uc-y)uc] = -VPc + V(77c(Vmc + (Vmc)^)) +—V- [/(rp)"^]; (1) 

Tp 

d,'^+{uc-V)'W = ■ (Vuc) + (Vuc)^■ (2) 

Tp 

Here, is the dynamic viscosity of the continuous solvent phase, Tfp the viscosity parameter for the EENE-P solute, 
Tp the polymer relaxation time, pc the solvent density. Pc the solvent bulk pressure, {VucY the transpose of 
the polymer-conformation tensor, I the identity tensor, f{rp) = {L^ — 3)/{L^ — rp) the EENE-P potential that ensures 
finite extensibility, rp = YTr{W) and L is the maximum possible extension of the polymers [26]. In the dispersed 
phase we just consider the NS equations 


Pd[dtUd + {ud-y)ud] = -VPrf + V(77rf(VMrf + (VMrf)^)) (3) 

with the shear viscosity of the dispersed phase. Immiscibility between the droplet phase and the matrix phase is 
introduced using the so-called “Shan-Chen” model [28, 29, 24] which ensures phase separation with the formation 
of stable interfaces between the two phases characterized by a positive surface tension a. As for the geometry used, 
we consider the simplest case where channels have a square cross-section H. The square cross-section is resolved 
with // X = 50 X 50 grid points, while the main channels and the side channel are resolved with 1150 and 250 grid 
points, respectively. We define, the fiow ratio as 2 = ^ = ^, the Capillary number as Ca = Reynolds number as 

Re = pcVcti !where and Qc = VcH^ are the fiow rates at the two inlets; vj and are the mean speeds in 

the inlets of the main and side channels, respectively, and cr is the interfacial tension. In all the numerical simulations, 

the viscous ratio between the disperse phase (r]j) and the matrix phase (rjc + rjp) is set to A = r]j /(rjc + rjp) = 1, with 

Rp/Vd = 0-4. The degree of viscoelasticity is computed from the Deborah number De = where Ni is the 

Ca 

first normal stress difference which develops in homogeneous steady shear [26, 30]. In the Oldroyd-B limit (L^ ^ 1), 
the Deborah number becomes which is clearly dependent on the ratio between the polymer relaxation 

time Tp and the emulsion time T^^ = • All the simulations described in this paper refer to a case with 

= 10^ (in LBM units). A proper tuning of the density gradients in contact with the wall allows for the modelling 
of the wetting properties which are chosen in such a way that the continuous phase completely wets the wall surfaces 
while the dispersed phase is non-wetting. 


RESULTS AND DISCUSSIONS 

In this section we report the results for the break-up of threads in a confined fiow-focusing device at fixed Capillary 
number. In a recent paper, Liu & Zhang [7, 8] performed LBM simulations and systematically characterized three 
typical flow patterns for different flow rate ratios at a fixed capillary number. At very low flow rate ratio 2, the droplets 
are formed at the cross-junction (DCJ) due to the squeezing mechanism [31]. Upon increasing Q, droplets are found to 
pinch-off downstream at the cross-junction (DC), forming a thread that becomes unstable after a distance from the 
cross-junction itself. Operatively, is computed as the distance between the break-up point and the corner up-stream 
at the cross-junction. As the flow rate ratio Q increases to a critical value, stable parallel flows (PE) are observed, where 
the three incoming streams co-flow in parallel to the downstream without pinching. In addition, the transitions from 
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(b) t = to 240Tem, Q = 2.5 
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(d) t = to 550Tem, Q = 2.5 


FIGURE 1. Dynamics and Break-up of liquid threads in flow-focusing geometry from LBM simulations. We report the time 
history for 4 representative time frames (a) prior to the first break-up, (b) in post-break-up conditions, (c) the droplet after the 
second break-up, and (d) finally the droplet prior to the fourth break-up. The case reported in this figure corresponds to Newtonian 
disperse phase (r]p = 0), Ca = 0.0056, Re = 0.015, X = I and Q = 2.5. The various snapshots show the time history of the liquid 
thread dynamics and break-up. We have chosen the emulsion time as unit of time (see text for details). 


DCJ to DC and from DC to PF are influenced by the Capillary number. As Ca increases, the threshold value of flow 
rate ratio at which the transition occurs decreases, and the width of the DC regime also decreases. Droplet formation 
and break-up in the DC regime is shown in Fig. 1 for Ca = 0.0056, Q = 2.5 for the case of a Newtonian Droplet in 
Newtonian matrix (r]p = 0 in Eq. (2)). Panels (a)-(d) show the droplet deformation prior to the first break-up at time 
t = to 150Tem, the droplet in post-break-up conditions at t = to X- 240 t^^, the droplet after the second break-up at 
r = ^0 + 410Tem, and finally the droplet prior to the fourth break-up at ^ + 550Tem, respectively. As we can see from 

the Fig. 1 the dispersed thread actually becomes unstable after a distance from the cross-junction and this distance 
increases as a function of time, which is another signature of the DC regime. 

To appreciate the effects of viscoelasticity, in Fig. 2 we show and compare the break-up process for different De. 
We report the density contours of the dispersed phase for three cases: (a) Newtonian matrix at time t = to 4 - 550Tem; (b) 
slightly viscoelastic matrix with De = 0.4 at time t = to 550Te„i; (c) viscoelastic matrix with Deborah number just 
above unity (De = 2.0) at time t = to 550t^^. Ah the other flow parameters are kept fixed, Ca = 0.0056, Re = 0.015, 
A = 1 and 2 = 2.5. We observe that as soon as the degree of viscoelasticity sensibly increases (De = 2.0), the breaking 
point of the liquid thread approaches the cross-junction (Panel (c)), whereas for De = 0.4 we only slightly perceive 
the effects of viscoelasticity. By comparing Panels (a) and (c) of Fig. 2 we also acknowledge a decrease in the droplet 
size, although such effect is not remarkable. By further increasing Q, also for the case of De = 2.0 the breaking 
point starts shifting downstream at the cross-junction. Ah these observations lend support to the fact that the effect of 
viscoelasticity results in a delayed transition from the DCJ regime to the DC regime. 


(a) t = to 550Tein, De^O 



(c) t = to 550Tem, De = 2.0 


FIGURE 2. Droplet formation in flow-focusing geometry at fixed Ca = 0.0056, Re = 0.015, A = 1 and Q = 2.5. Three cases are 
compared: (a) Newtonian case (De = 0.0) at time t = to-\- 550Tein; (b) slightly viscoelastic case with De = 0.4 at time t = to-\- 550Tem; 
(c) viscoelastic case with De = 2.0 at time t = to + 550Tem. 




CONCLUSIONS 

Microfluidic technologies offer the possibility to generate small fluid volumes of dispersed phases (droplets) in 
continuous phases. One of the most common droplet generator is the flow-focusing device, based on the focusing 
a stream of one liquid in another co-flowing immiscible liquid. The conflnement that naturally accompanies these 
devices has an impact on droplet deformation and break-up, which are signiflcantly different from those of unbounded 
droplets. Moreover, in real processing conditions, relevant constituents have commonly a viscoelastic -rather than 
Newtonian- nature. In this paper we have presented results based on numerical simulations with the "lattice Boltzmann 
models" (LBM) to highlight the non trivial role played by conflnement and non-Newtonian effects. We have worked 


in the operating conditions previously studied by other authors [7, 8, 32], where Newtonian droplets in a Newtonian 
matrix are formed either at the cross-junction, due to the squeezing mechanism [31], or downstream at the cross¬ 
junction (DC), forming a thread that becomes unstable after a distance from the cross-junction. The transition 
between these two regimes has been found to be affected by matrix viscoelasticity. Future investigations on the role 
of droplet viscoelasticity are definitively worth of being considered. We also remark that the results presented in this 
paper are obtained for a fixed L^, i.e. for a fixed maximum elongation of the polymers. The parameter is actually 
related to the extensional viscosity of the polymers [26, 30]. In the spirit of the work that we recently developed 
for characterizing droplet deformation and break-up in confined shear flow [25], it would be of extreme interest to 
study the impact of a change in the finite extensibility parameter for the geometries studied in this paper [13, 16]. 
Finally, we wish to underscore the role played by numerical simulations: in real experiments all the various processes 
(hydrodynamic forces, viscoelastic effects, confinement effects) occur at the same time and it is next to impossible to 
separately quantify their relative importance. The opportunity offered by numerical studies is to allow for a systematic 
analysis of each of the effects separately. 
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